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Abstract. In this paper, I test by numerical simulations the results of Del Popolo & 

m . 

■ Gambera (1998), dealing with the extension of Chandrasekhar and von Neumann's anal- 
ysis of the statistics of the gravitational field to systems in which particles (e.g. stars, 

^ | galaxies) are inhomogeneously distributed. The paper is an extension of that of Ahmad 
,H \ & Cohen (1974), in which the authors tested some results of the stochastic theory of 

■ dynamical friction developed by Chandrasekhar & von Neumann (1943) in the case of 
^ ■ homogeneous gravitational systems. It is also a continuation of the work developed in Del 

Popolo (1996a, b), which extended the results of Ahmad & Cohen (1973), (dealing with 
5t i the study of the probability distribution of the stochastic force in homogeneous gravi- 
tational systems) to inhomogeneous gravitational systems. Similarly to what was done 
by Ahmad & Cohen (1974) in the case of homogeneous systems, I test, by means of the 
evolution of an inhomogeneous system of particles, that the theoretical rate of force fluc- 
tuation dF/dt describes correctly the experimental one, I find that the stochastic force 
distribution obtained for the evolved system is in good agreement with the Del Popolo 
& Gambera (1998) theory. Moreover, in an inhomogeneous background the friction force 
is actually enhanced relative to the homogeneous case. 
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1. Introduction 

The study of the statistics of the fluctuating gravitational force in infinite homogeneous 
systems was pioneered by Chandrasekhar & Von Neumann in two classical papers (Chan- 
drasekhar & Von Neumann 1942, 1943 hereafter CN43) and in several other papers by 
Chandrasekhar (1941, 1943a, 1943b, 1943c, 1943d, 1943e, 1944a and 1944b). The analysis 
of the fluctuating gravitational field, developed by the authors, was formulated by means 
of a statistical treatment. In their papers Chandrasekhar & Von Neumann considered a 
system in which the stars are distributed according to a uniform probability density, no 
correlation among the positions of the stars is present and where the number of stars 
constituting the system tends to infinity while keeping the density constant. Two distri- 
butions are fundamental for the description of the fluctuating gravitational field: 

1) W(F) which gives the probability that a test star is subject to a force F in the range 
F, F+ dF; 

2) W(F, f ) which gives the joint probability that the star experiences a force F and a 
rate of change f, where f = dF/dt. 

The first distribution, known as Holtsmark's law (Holtsmark 1919), in the case of a 
homogeneous distribution of the stars, gives information only on the number of stars 
experiencing a given force but it does not describe some fundamental features of the fluc- 
tuations in the gravitational field such as the speed of the fluctuations and the dynamical 
friction. These features can be described using the second distribution W(F,f). Hence, 
for the definition of the speed of fluctuations and of the dynamical friction one must 
determine the distribution W(F, f). Information on dynamical friction can be obtained 
from the moments of this last distribution. As stressed by Chandrasekhar & von Neu- 
mann (1943), for a test star moving with velocity v in a sea of field stars characterized 
by a random probability distribution of the velocities, $(u), we may write: 



where V represents the velocity of a typical field star relative to the one under consid- 
eration, u denotes the velocity of a field star. This asymmetry of the distribution of the 
relative velocities produces, as shown by CN43, a deceleration of the test star in the 
direction of motion. This effect is known, from the Chandrasekhar papers, as "dynamical 
friction^ . As quoted above, some information on dynamical friction can be obtained by 
means of the first moment of f. As shown by CN43: 



(V) = (u) - v = -v 



(1) 
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where m is the mass of a field star, ni is the local density, and B{f3) is definited as: 

. , In H((3)d(3 

B ^= 3 V^r- 1 (4) 

where 

roo 3 

H{(3) = — e~ifr xsin(x)dx (5) 
Jo 

and p = \F\/Qh where Qh = (jk) 2 ^ ZftGMn 2 / 3 . These equations show that the amount 
of acceleration in the direction of — v when v • F < is greater than that in the direction 
+v, when v • F > 0: then the star suffers a deceleration, the a priori probability that 
v • F > being equal to the probability that v • F < 0. 

Several authors have stressed the importance of stochastic forces and in particular dynam- 
ical friction in determining the observed properties of clusters of galaxies (White 1976; 
Kashlinsky 1986,1987) while others studied the role of dynamical friction in the orbital de- 
cay of a satellite moving around a galaxy or in the merging scenario (Bontekoe & van Al- 
bada 1987; Scguin & Dupraz 1996; Dominguez-Tenreiro & Gomez-Flechoso 1998). Chan- 
drasekhar's theory (and in particular his classical formula (see Chandrasekhar 1943b)) 
is widely employed to quantify dynamical friction in a variety of situations, even if the 
theory developed is based on the hypothesis that the stars are distributed uniformly and 
it is well known that in stellar systems, the stars are not uniformly distributed, (Elson 
et al. 1987; Wybo & Dcjonghe 1996; Zwart et al. 1997) in galactic systems as well, the 
galaxies arc not uniformly distributed (Peebles 1980; Bahcall & Soncira 1983; Sarazin 
1988; Liddle, & Lyth 1993; White et al. 1993; Strauss & Willick 1995). It is evident that 
an analysis of dynamical friction taking account of the inhomogeneity of astronomical 
systems can provide a more realistic representation of the evolution of these systems. 
Moreover from a pure theoretical standpoint we expect that inhomogeneity affects all 
the aspects of the fluctuating gravitational field (Antonuccio & Colafrancesco 1994; Del 
Popolo & Gambera 1996, 1997; Del Popolo et al. 1996; Gambera 1997). Firstly, the Holts- 
mark distribution is no longer correct for inhomogeneous systems. For these systems, as 
shown by Kandrup (1980a, 1980b, 1983), the Holtsmark distribution must be substituted 
with a generalized form of the Holtsmark distribution characterized by a shift of W(F) 
towards larger forces when inhomogeneity increases. This result was already suggested 
by the numerical simulations of Ahmad & Cohen (1973, 1974). Hence when the inhomo- 
geneity increases the probability that a test particle experiences a large force increases, 
secondly, W(F,f) is changed by inhomogeneity. Consequently, the values of the mean 
life of a state, the first moment of f and the dynamical friction force are changed by 
inhomogeneity with respect to those of homogeneous systems. 

This paper must be intended as the continuation of Del Popolo & Gambera (1998) paper, 
in which the study of the effects of inhomogeneity on the distribution functions of the 
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stochastic forces and on dynamical friction was performed. 

As anticipated in Del Popolo & Gambera (1998), the next task to perform was to test 
the result of the Del Popolo & Gambera (1998) paper against N-body simulations, which 
is the object of the present paper. The third step (to be developed in a future paper) 
should be that of finding a formula that describes dynamical friction in homogeneous 
and inhomogeneous systems only on the basis of statistical theory. 

Before continuing we want to stress that when we speak of inhomogeneity we refer to 
inhomogeneity in position distribution and not to that of velocity distribution. Our work 
follows the spirit of Kandrup's (1980) in the sense that we are interested in the effect of 
a non-uniform distribution in the position of stars on the distributions of the stochastic 
force. 

The plan of the paper is the following: in Sec. 2 we review the calculations and 
formulas needed to obtain (^-)f,v m the case of inhomogeneous systems. In Sec. 3, 1 show 
how numerical experiments are performed and they are compared with the theoretical 
results of Del Popolo & Gambera (1998). Finally, in Sec. 4, I draw my conclusions. 



2. Force derivative and dynamical friction 

The introduction of the notion of dynamical friction is due to CN43. In the stochastic 
formalism developed by CN43 the dynamical friction is discussed in terms of f : 

3F-v 



27T GmnB([3) 



(6) 



3 

where 

As shown by CN43, the origin of dynamical friction is due to the asymmetry in 
the distribution of relative velocities. If a test star moves with velocity v in a spherical 
distribution of field stars, namely </>(u), then we have that: 

V = u=v = -v (8) 

The asymmetry in the distribution of relative velocities is conserved in the final Eq. 10 • 
In fact from Eq. © we have: 

(CN43). This means that when v F > then ^ > 0, while when v • F < then 
-^r < 0. As a consequence, when F has a positive component in the direction of v, 
|F| increases on average; while if F has a negative component in the direction of v, |F| 
decreases on average. Moreover, the star suffers a greater amount of acceleration in the 
direction — v when v • F < than in the direction +v when v • F > 0. 
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In other words the test star suffers, statistically, an equal number of accelerating and 
decelerating impulses. The modulus of deceleration being larger than that of acceleration 
the star slows down. 



Chandrasekhar & von Neumann's analysis was extended to inhomogeneous systems 
in Del Popolo & Gambera (1998), where W(F, f ) and its first momentum were calculated. 
Supposing that the distribution function r is given by: 



-^0' 2 (M)|u| : 



(10) 



where a is a constant that can be obtained from the normalization condition for r, j 
a parameter (of dimensions of velocity -1 ), ip an arbitrary function, u the velocity of a 
field star. In other words we assumed, according to CN43 and Chandrasekhar & von 
Newmann (1942), that the distribution of velocities is spherical, i.e. the distribution 
function is ip(u) = ^(^(Af^ul 2 ), but differently from the quoted papers we supposed 
that the positions are not equally likely for stars, that is, the stars are inhomogeneously 
distributed in space. A lengthy calculation led us (see Del Popolo & Gambera 1998 for a 
derivation and the meaning of symbols) to find that the first moment of f is given by: 



f = 

where: 
A(p) = 

B(p) = 



A(p) ■ B(p) 3 -p 



a~GML((3) 
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The results obtained by us for an inhomogeneous system are different [see Eq. Hilt ]. a s 
expected, from those obtained by CN43 for a homogeneous system (CN43 - Eq. 105 or Eq. 
(©). At the same time it is very interesting to note that for p = (homogeneous system) 
our result coincides, as is obvious, with the results obtained by CN43. In inhomogeneous 
systems, Eq. (0 can be written, using Eq. ({11(1 . as: 



,d|F|. v a u~GML{B) vF 



(16) 
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In order to check the validity of the quoted relation (Eq. <|16l) h I have performed numeri- 
cal experiments. This was done similarly to Del Popolo (1996b) by evolving (now) 100000 
points (stars) acting under their mutual gravitational attraction. From the evolved posi- 
tions and velocities of the stars, (^r) was computed as a function of velocity and force, 
similarly to Ahmad & Cohen (1974), and then compared with Eq. (|16fl as I shall describe 
in the following. 1 

3. N-body experiments 

To calculate the stochastic force in an inhomogeneous system, I used an initial config- 
uration in which particles were distributed according to a truncated power-law density 
profile: 

p(r) = p (^) P 0<r<R (17) 

(see Kandrup 1980, Del Popolo 1996a). If the velocity distribution is everywhere isotropic 
then the equation relating the configuration space density p(r) to the phase space density 
f(E) is: 



p(r) = 4tt / ^2 [E - U{r)]f(E)dE (18) 

JU(r) 

where U{r) is the potential (normalized to zero at infinity). Eq. I|18|) may be converted 
into an Abel integral equation and inverted, giving the phase space density: 

(Eddington 1916; Binney & Tremaine 1987). The initial conditions were generated from 
the distribution function that can be obtained from Eq. (|19|) assuming a cut-off radius 
R = 1, the mass of the system M=l,r = 0.15 2 and G — 1. All the particles had equal 
mass. To have a system whose total mass is contained in a unitary sphere, Eq. i|18|) was 
renormalized and consequently also the potential of the system which is obtained from 
Eq. (|18(l through Poisson's equation. The system of 100000 particles was evolved over 
150 dynamical times using a tree N-body code (Hernquist 1987). During the evolution 
of the system the essential quantities such as position, force, etc., for each test point of 
the system was sampled every ^ of a dynamical time (see Del Popolo 1996b, Ahamad 
& Cohen 1973 for details). The average of the is a function of velocity, v, force, 

1 In a inhomogeneous system, similarly way to what happens in a homogeneous system, f 
depends on v, F and 9 (the angle between v and F) while unlike homogeneous systems, f is a 
function of the inhomogeneity parameter p. The dependence of f on p is not only due to the 
functions A(p), B(p) and to the density parameter a but also to the parameter /3 = |F|/Qh. 
In fact in inhomogeneous systems the normal field Qh is given by Qh = GM(aB(p)/2) 2 ^ 3 ~ p \ 

clearly dependent on p. 

2 This is the value I used, remember however that the distribution is scale-free 
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Fig. 1. The average value of the time rate of change of the magnitude of the force as the 
function of the velocity. The solid line refers to the homogeneous case (Chandrasekhar & 
von Neumann 1943). The dotted and dashed line refers to the cases p — 2.5 and p = 4, 
respectively (see Eq. 26). Crosses represent the experimental points. 

F, and the angle between them. The test of Eq. (|XfSf) was performed in a similar way 
to that of Ahmad & Cohen (1974), namely by integrating out two of the variables and 
examining (^jr^) against the remaining one (see Ahmad & Cohen 1974). As in Ahmad & 
Cohen (1974), (^^)ipi indicates (^^) after integrating out the angle and velocity, while 




Fig. 2. The average value of the time rate of change of the magnitude of the force as the 
function of the force. The solid line refers to the homogeneous case (Chandrasekhar & 
von Neumann 1943). The dotted and dashed line refers to the cases p = 2.5 and p = 4, 
respectively (see Eq. 19). Crosses represent the experimental points. 

d I F I d I F I 

( dt )|v| i s (~3t ) a ft cr integrating out the angle and force. In integrating out the angle, 
one cannot average over the entire range to 7r, since that would give zero. Instead, the 
cosine was averaged from to 7r/2, and in the numerical experiments only those particles 
having a cosine in the quoted range were used. In order not to waste the statistics for 
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half the particles the same trick of Ahmad & Cohen (1974) was used, namely when the 
cosine is in the range of 7r/2 to n, the sign of is changed and it is counted in the same 
statistics. This corrensponds to assume that the cosine between v and F is uniformly 
distributed, which is what is found in numerical experiments. The average value of the 
cosine between and ir/2 is 1/2. 

For a general distribution (|v|) can be calculated as usual: 

and (^gp-) can be written, in units of — a^p GM (|v|) , as: 

In the particular case of a Maxwellian distribution for velocities: 

= J^e X p(~] 2 \v\ 2 ) (22) 
where j 2 = 2< ^ > , so that: 

<|v|) = -i: (23) 

7T ' J 

we have that: 



,d|F|. 2— _b_ a-GML(P) 
— ■ A(p) ■ B(p) s-f 

2 J ffW 



that expressed in units of ^; 3 . a 3 -pGM , then Eq. (|2*4"|) becomes: 

<ffi )|p| .«*.A (25) 

Similarly to Ahamd & Cohen (1974), since to integrate out the force from Eq. (|16|) 
one has a divergent result, I consider only particles up to a certain maximum value of 
the force, /3 max : for example in the case p = /3 max = 10.7, that involves 97% of the 
particles. As observed by Ahmad & Cohen (1974), any cutoff of the force can be used 
as long as it is taken into account in both the experiment and the analytic evaluation of 
(^)|v|- Defining: 

Bi{P) = L(P \- P (26) 
and 

we find, in units of (3)^ (Bl(/3))a^GM, that: 

(^)|v|=^(p)-B(p)^-« (28) 
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The results of calculation and numerical experiments are plotted in Fig. 1 and Fig. 
2. In Fig. 1, I plot the average value of the time rate of change of the magnitude of 
the force as a function of the velocity. The solid line refers to the homogeneous case 
while the dotted and dashed lines refer to the cases p = 2.5 and p — 4, respectively. 
Crosses represent the experimental points. As shown, experimental points follow a linear 
relationship and there is a good agreement with the theoretical prediction, (Eg. I21f> . In 
Fig. 2, 1 plot the average value of the time rate of change of the magnitude of the force as 
the function of the force. As in the previous figure, the solid line refers to the homogeneous 
case while the dotted and dashed line refers to the cases p = 2.5 and p = 4, respectively. 
In this case, the dependence is no longer linear: it behaves like B((3) in the homogeneous 
case. A comparison with numerical experiments shows that there is a good agreement 
with the theoretical prediction, (Eq. I28[l . The situation described in Ahmad & Cohen 
(1974), that the experimental data were somewhere in between the theoretical curves for 
the one-particle and the infinite-particle case, is no longer present and the agreement 
is better, now. This is due to the larger number of particles used in the simulations. 
The plots show that Chandrasekhar & von Neumann's theory of dynamical friction in 
gravitational systems gives a good description of experimental data (solid line and data), 
and so does the generalization of the quoted theory to inhomogeneous systems (dotted 
line, dashed line, and respective data). In inhomogeneous systems, Chandrasekhar's result 
which relates the frictional force only to the local properties of the background at the 
position of the object, is no longer true, and friction depends on the global structure of the 
system. This point is in agreement with Maoz (1993), who showed that in inhomogeneous 
media the friction, unlike Chandrasekhar's formula, depends on the global structure of 
the entire mass density field. 



4. Conclusions 

In this paper, I tested by numerical simulations the results of the Del Popolo & Gambera 
(1998) paper, dealing with the average value of the time rate of change of the magnitude 
of the stochastic force in inhomogeneous gravitational systems. In agreement with Ahmad 
& Cohen (1974), the stochastic theory of dynamical friction developed by Chandrasekhar 
& von Neumann (1943), in the case of homogeneous gravitational systems, gives a good 
description of the results of numerical experiments. The stochastic force distribution ob- 
tained for inhomogeneous systems, obtained by Del Popolo & Gambera (1998), is also in 
good agreement with the results of numerical experiments. Finally, in an inhomogeneous 
background the friction force is actually enhanced relative to the homogeneous case. 
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